%==========================================================================
%=  This file is part of the Sound Restoration Project
%=  (c) Copyright Industrial Mathematics Institute
%=                University of South Carolina, Department of Mathematics
%=  ALL RIGHTS RESERVED
%=
%=  Author: Borislav Karaivanov
%==========================================================================

%==========================================================================
% List of the files on which this procedure depends:
%
% curve fitting toolbox
%
%==========================================================================

%==========================================================================
% The function "visualizeSoundCalibrationSingleMask" visualizes the
% elements of the single-mask sound calibration process (averages,
% calibration mask, calibrated signals).
% INPUT: "origSignalArr" is a column vector holding the signal to be used
% in the visualization of the calibration process.
% "aveArr" is a column vector holding the average of the signals
% contributing to the calibration mask.
% "sweeps" is a non-negative integer specifying the number of smoothing
% passes to be applied to the raw average of the given signals.
% "span" is a positive odd integer specifying the span of the stencil to be
% used for moving average smoothing of the average signal.
% "degOfPolyFit" is an optional, non-negative integer specifying the degree
% of the polynomial to be fit to the average signal. This polynomial is to
% be sampled and used as a calibration mask. If this parameter is
% specified, then the values of the previous two parameters are ignored.
% OUTPUT: None.
%==========================================================================
function visualizeSoundCalibrationSingleMask(origSignalArr, aveArr, ...
    sweeps, span, degOfPolyFit)

% Use the existing global variables to determine whether visualization is
% required, and what common prefix and file extension are to be used for
% the names of the image files saving snapshots of figures.
global globalDoVisualize globalShotFileNamePrefix globalShotFileExt
if (globalDoVisualize ~= true)
    return;
else
    shotFileNamePrefix = globalShotFileNamePrefix;
    shotFileExt = globalShotFileExt;
end

% If any input argument is not specified (i.e., is empty) then assign it a
% default value so that the visualization can be made.
if (isempty(degOfPolyFit) == true)
    degOfPolyFit = 2;
end
if (isempty(sweeps) == true)
    sweeps = 1;
end
if (isempty(span) == true)
    span = 15;
end

% Allocate memory for the averages.
signalLength = size(aveArr, 1);
aveAllArr = zeros(signalLength, sweeps + 1);

% Initialize with the average of all signals.
aveAllArr(:, 1) = aveArr;

% Compute the desired number of smoothed averages using moving average.
for k = 1:sweeps
    aveAllArr(:, k + 1) = smooth(aveAllArr(:, k), span, 'moving');
end

% Compose the current common prefix.
commonPrefix = [shotFileNamePrefix '_Sweeps' num2str(0) '_Span' num2str(span)];

% Plot the raw average.
figHandle = figure;
plot(aveAllArr(:, 1));
xlabel('Average signal (never smoothed)');
axis tight;
% Maximize the figure to cover the whole screen.
set(figHandle, 'Position', get(0, 'ScreenSize'));
% Save a snapshot.
shotFileName = [commonPrefix '_Ave.' shotFileExt];
print(gcf, ['-d' shotFileExt], shotFileName);
% Close the figure.
close(figHandle);

% Plot a signal and its calibrated version.
figHandle = figure;
plot((1:signalLength), origSignalArr, 'b', (1:signalLength), origSignalArr./aveAllArr(:, 1), 'r');
xlabel('Original signal in blue, calibrated signal in red');
title('Calibration mask based on the average signal (never smoothed)');
axis tight;
% Maximize the figure to cover the whole screen.
set(figHandle, 'Position', get(0, 'ScreenSize'));
% Save a snapshot.
shotFileName = [commonPrefix '_OrigAndCalib.' shotFileExt];
print(gcf, ['-d' shotFileExt], shotFileName);
% Close the figure.
close(figHandle);

if (sweeps > 0)
    % Compose the current common prefix.
    commonPrefix = [shotFileNamePrefix '_Smooth' num2str(1) '_Span' num2str(span)];
    
    % Plot the first smoothed average.
    figHandle = figure;
    plot(aveAllArr(:, 2));
    xlabel(['Average signal (smoothed once with moving average filter of span ' num2str(span) ')']);
    axis tight;
    % Maximize the figure to cover the whole screen.
    set(figHandle, 'Position', get(0, 'ScreenSize'));
    % Save a snapshot.
    shotFileName = [commonPrefix '_Ave.' shotFileExt];
    print(gcf, ['-d' shotFileExt], shotFileName);
    % Close the figure.
    close(figHandle);
    
    % Plot a signal and its calibrated version.
    figHandle = figure;
    plot((1:signalLength), origSignalArr, 'b', (1:signalLength), origSignalArr./aveAllArr(:, 2), 'r');
    axis tight;
    xlabel('Original signal in blue, calibrated signal in red');
    title(['Calibration mask based on one pass of moving average with span of ' num2str(span)]);
    % Maximize the figure to cover the whole screen.
    set(figHandle, 'Position', get(0, 'ScreenSize'));
    % Save a snapshot.
    shotFileName = [commonPrefix '_OrigAndCalib.' shotFileExt];
    print(gcf, ['-d' shotFileExt], shotFileName);
    % Close the figure.
    close(figHandle);
end

for k = 2:sweeps
    % Compose the current common prefix.
    commonPrefix = [shotFileNamePrefix '_Smooth' num2str(k) '_Span' num2str(span)];
    
    % Plot the rest of the smoothed averages.
    figHandle = figure;
    plot(aveAllArr(:, k + 1));
    xlabel(['Average signal (smoothed ' num2str(k) ' times with moving average filter of span ' num2str(span) ')']);
    axis tight;
    % Maximize the figure to cover the whole screen.
    set(figHandle, 'Position', get(0, 'ScreenSize'));
    % Save a snapshot.
    shotFileName = [commonPrefix '_Ave.' shotFileExt];
    print(gcf, ['-d' shotFileExt], shotFileName);
    % Close the figure.
    close(figHandle);
    
    % Plot a signal and the rest of its calibrated versions.
    figHandle = figure;
    plot((1:signalLength), origSignalArr, 'b', (1:signalLength), origSignalArr./aveAllArr(:, k + 1), 'r');
    xlabel('Original signal in blue, calibrated signal in red');
    title(['Calibration mask based on ' num2str(k) ' passes of moving average with span of ' num2str(span)]);
    axis tight;
    % Maximize the figure to cover the whole screen.
    set(figHandle, 'Position', get(0, 'ScreenSize'));
    % Save a snapshot.
    shotFileName = [commonPrefix '_OrigAndCalib.' shotFileExt];
    print(gcf, ['-d' shotFileExt], shotFileName);
    % Close the figure.
    close(figHandle);
end


% Compute a calibration mask by a polynomial fit of desired degree by
% fitting the average signal with a polynomial of the desired degree, and
% sampling the found polynomial to produce the calibration mask.
xArr = (1:signalLength).';
curveFunct = fit(xArr, aveAllArr(:, 1), ['poly' num2str(degOfPolyFit)]);
fitCalibMaskArr = curveFunct(xArr);

% Plot the smoothed and fitted masks.
figHandle = figure;
plot(xArr, aveAllArr(:, 1), '-b');
hold on;
plot(xArr, fitCalibMaskArr, '-k', xArr, aveAllArr(:, end), '-r', 'LineWidth', 1);
xlabel(['Average in blue, smoothed average in red, polynomial fit of degree ' num2str(degOfPolyFit) ' in black']);
axis tight;
% Maximize the figure to cover the whole screen.
set(figHandle, 'Position', get(0, 'ScreenSize'));
hold off;
% Save a snapshot.
shotFileName = [commonPrefix '_Deg' num2str(degOfPolyFit) '_Masks.' shotFileExt];
print(gcf, ['-d' shotFileExt], shotFileName);
% Close the figure.
close(figHandle);

% Plot a signal and its calibrated version.
figHandle = figure;
plot((1:signalLength), origSignalArr, 'b', (1:signalLength), origSignalArr./fitCalibMaskArr, 'r');
axis tight;
xlabel('Original signal in blue, calibrated signal in red');
title(['Calibration mask based on polynomial fit of degree ' num2str(degOfPolyFit)]);
% Maximize the figure to cover the whole screen.
set(figHandle, 'Position', get(0, 'ScreenSize'));
% Save a snapshot.
shotFileName = [shotFileNamePrefix '_PolyDeg' num2str(degOfPolyFit) '_OrigAndCalib.' shotFileExt];
print(gcf, ['-d' shotFileExt], shotFileName);
% Close the figure.
close(figHandle);

% Plot the difference between the fitted calibration mask and the average
% signal.
figHandle = figure;
plot(xArr, fitCalibMaskArr - aveAllArr(:, 1), '-b', xArr, aveAllArr(:, end) - aveAllArr(:, 1), '-r');
xlabel(['Differences between: polynomial fit of degree ' num2str(degOfPolyFit) ' and average in blue, smoothed average and average in red']);
axis tight;
% Maximize the figure to cover the whole screen.
set(figHandle, 'Position', get(0, 'ScreenSize'));
% Save a snapshot.
shotFileName = [commonPrefix '_Deg' num2str(degOfPolyFit) '_MaskDiff.' shotFileExt];
print(gcf, ['-d' shotFileExt], shotFileName);
% Close the figure.
close(figHandle);

return;
% end of the function "visualizeSoundCalibrationSingleMask"
